%% peformance variables
% 1) Fs = 2000;
% 2) Usage of several trials per movement and not just one (currently)
% 3) Small shift since the trigger event
% 3) Selected frequencies
% 4) Potentially NCV

%% segement all the files

FileList = {
   
   'mohEEG111.mat'
  % 'Marion12_2504.mat'
%     'Tim12_2804.mat'
    };

%FileList = { 'session2.acq.mat' };
%Fs = 2000; default
Fs = 2000;

%window = Fs;
window = 6000; % 3 seconds
X = []; % contains the trials from all classes
Y = []; % contains the classification from all classes

frequencies = {[20 100] [280 320] [580 620]};
%frequencies = {[20 80] [80 150] [150 250] [250 300] [250 300] [300 350]  [350 500] [500 700] };

Nfreq = length(frequencies);

for curFileIx = 1:length(FileList)
    %load(['D:/dropbox/Gipsa-work/data_EMG/franck/' FileList{curFileIx}]);
    %load(['D:\Dropbox\Gipsa-work\data_EMG\delmas\' FileList{curFileIx}]);
    load(['D:\temp\subjects\' FileList{curFileIx}]);
    disp(['Subject file:' FileList{curFileIx}]);
    
    datar = data(:,1:size(data,2)-2); % only data
    
    %trial_channel = data(:,8);
    trial_channel = data(:, size(data,2)-1);
    
    %class_separation_channel = data(:,9);
    class_separation_channel = data(:, size(data,2));
    
    classN = length(count_sequence(class_separation_channel,5));
    disp(['Classes/Movements detected: ' int2str(classN)]);
    
    
    trialsN = length(count_sequence(trial_channel,5));
    disp(['Trials per class detected: ' int2str(trialsN)]);
    
    
    for m = 1:classN
        
        %classes should be supplied to this function
        %select only the trials for this class
        %"best_index" - in 1 shows the postion of a trial, "best_trigger" -
        %again shows the postion but in the range 1 .. n
        [best_index,best_trigger] = generate_trials_postitions_using_discrete_channel(m, trial_channel,class_separation_channel);

        %plot_trials(data,FileList{curFileIx},best_trigger);

        % cut signal
        trials_per_class = zeros(size(datar,2),window,length(best_index),Nfreq);
        disp(['Trials: ' int2str(length(best_index)) ' class ' int2str(m)]);
        
        for f=1:Nfreq

            [b,a] = build_filter('Fs',Fs,'HP',frequencies{f}(1),'LP',frequencies{f}(2));
            current_filtered_data = filtfilt(b,a,datar);

            for t=1:length(best_index)
                %calculate a range around the trigger
                shift = 0;% 100 with 2000 sampling frequency this is 50 ms

                selected_interval = best_index(t)+shift:best_index(t)+shift+window-1;
                
                trials_per_class(:,:,t,f) = current_filtered_data(selected_interval,:)';      
            end;
        end;

        X = cat(3,X,trials_per_class); % put all trials together - from all classes/files/movements
        
        Y = cat(1,Y,ones(length(best_index),1) * m);
    end;
    
    disp('--------------------------------');
    
end

%% Classification over X and Y

%shuffuling
NTrial = size(X,3); % trials equals movements
disp(['Total trials from all classes and all subjects: ' int2str(NTrial)]);
ix = randperm(NTrial);

% permutate the data for classification
X = X(:,:,ix,:);
Ys = Y(ix); % it takes Y in the permutated form

% cross validation
NCV = 10;
% MDM
disp('---------- MDM --------------');
out = cross_valid(X,Ys,NCV,@mdm_freq); % COV is produced from X and Ys is produced from Y
disp('Confusion matrix :');
disp(confusionmat(Ys,out));
disp(['Accuracy : ' num2str(100*mean(Ys==out)) ' %']);
